---
title: "Fairweather Cosmo: Supplemental Info"
output: html_notebook
---


```{r setup, warning=F, message=F}

# Load pacman
library(pacman)

# Use pacman to load libraries
p_load(tidyverse, haven, lmerTest, xtable, janitor, knitr, huxtable, marginaleffects,
       ggplot2, ggpubr, ggsci, countrycode, readxl, texreg, dotwhisker, ggrepel,
       here, psych, ggcorrplot, Hmisc, xtable, conflicted)

# Conflicts 
conflicts_prefer(psych::alpha)
conflicts_prefer(dplyr::filter)
conflicts_prefer(dplyr::summarize)
conflicts_prefer(lme4::lmer)

# Read the WVS data
wvs_latam <- read_csv("wvs_latam.csv")

# Make the var names lowercase
wvs_latam <- clean_names(wvs_latam)



```



## Correlations 

Corrs between main variables:

```{r corr}

# variables 
corr_vars <- wvs_latam |> 
  select(imm_index, natpride, close_world, close_world_z, close_country,
         close_country_z)

# correlations
corrs <- corr_vars |> 
  na.omit() |> 
  cor(method = "spearman") |> 
  round(digits = 2)

# labels
colnames(corrs) <- c("Immigration Index", "National Pride",
               "Cosmopolitianism \n(Close to World)", "Cosmo Z-Score", 
               "Close to Country", "Close to Country Z-Score")

rownames(corrs) <- c("Immigration Index", "National Pride",
               "Cosmopolitianism", "Cosmo Z-Score", 
               "Close to Country", "Close to Country Z-Score")

# plot
ggcorrplot(corrs,
           lab = T, 
           type = "lower")

# save
ggsave(file = "figures/corrs1.pdf", height = 4, width = 6)

```

National pride and other out-group variables:

```{r corr-2}

# select vars 
vars2 <- wvs_latam |>
  select(
    natpride, close_country, trust_other_nat, diff_race_neigh, diff_lang_neigh, 
    imm_diversity, imm_socconflict, imm_terrorism, imm_unemployment, imm_development
    ) |> 
  data.frame()

# labels
corr_labs <- c("National Pride", "Close Country", "Trust Other Nationality", "Neighbor Different Race", "Neighbor Different Language", "Imm. Strengthens Diversity", "Imm. Reduces Conflict", "Imm. Decreases Terrorism", "Imm. Decreases Unemployment", " Imm. Good for Development")


# correlation
corrs2 <- vars2 |> 
  na.omit() |>
  cor(method = "spearman")

# p values 
corr2_p <- corrs2 |> 
  na.omit() |> 
  cor_pmat(method = "spearman")

# label matrices
rownames(corrs2) <- corr_labs
colnames(corrs2) <- corr_labs
rownames(corr2_p) <- corr_labs
colnames(corr2_p) <- corr_labs

# plot
ggcorrplot(corrs2, 
           lab = T,
           type = "lower") 
           #p.mat = corr2_p 

# save
ggsave(file = "figures/corrs2.pdf", height = 6, width = 8) 


```


# Recoding/Merging for Models 

```{r shaping}
# Make Education a Factor
wvs_latam <- wvs_latam %>% 
  mutate(edu_cat = case_match(edu,
                              1:2 ~ "Primary", 
                              1:4 ~ "Secondary",
                              5:8 ~ "Post-Secondary"))

wvs_latam$edu_cat <- factor(wvs_latam$edu_cat, levels = c("Primary", 
                                              "Secondary", 
                                              "Post-Secondary"))


# Log the Income Level Variable 
wvs_latam$income_log <- log(wvs_latam$income_lvl)

# Log Age 
wvs_latam$ln_age <- log(wvs_latam$age)

# Make a couple of variables factors
wvs_latam$unemp <- factor(wvs_latam$unemp) #unemployment
wvs_latam$female <- factor(wvs_latam$female) #female
wvs_latam$urban <-  factor(wvs_latam$urban) # urban rural


```


# Models With Z Scores

```{r z score mods}

# Main Effect of Cosmo
m1z <- lmer(imm_index ~ close_world_z + edu_cat + trust_generalized + income_log + urban +
             ln_age + female + unemp + refugee_change + ln_gdppc +  (1| iso3c), data = wvs_latam) 



# Conditional Effect of Cosmopolitanism 
m2z <- lmer(imm_index ~ close_world_z * refugee_change + edu_cat + trust_generalized + 
             income_log + urban + unemp + ln_age +  female + unemp + ln_gdppc +  (1| iso3c), data = wvs_latam) 


# Main Effect of National Pride
m3z <- lmer(imm_index ~ close_country + edu_cat + trust_generalized + income_log + ln_age + female + 
             unemp + urban + refugee_change + ln_gdppc + (1| iso3c), data = wvs_latam) 


# Conditional Effect of National Pride
m4z <- lmer(imm_index ~  close_country_z * refugee_change +  edu_cat + trust_generalized + income_log + urban + ln_age + female + unemp + ln_gdppc + female + (1| iso3c), data = wvs_latam) 


# Coefficient Labels 
coef_labs_z <- list(
                  "close_world_z" = "Cosmopolitanism Z-Score", 
                   "close_world_z:refugee_change" = "Cosmopolitanism Z-Score $\\times$ Refugee Change",
                  "close_country_z" = "Close Country Z-Score",
                  "close_country_z:refugee_change" = "Close Country Z-Score $\\times$ Refugee Change",
                  "trust_generalized" = "Interpersonal Trust",
                  "edu_catSecondary" = "Secondary Edu", 
                  "edu_catPost-Secondary" = "Post-Secondary Edu", 
                  "ln_age" = "Age (logged)",
                  "urban1" = "Urban",
                  "unemp1" = "Unemployed",
                  "female1" = "Female",
                  "income_log" = "Individual Income", 
                  "refugee_change" = "Refugee Change", 
                   "ln_gdppc" = "Log GDP/Capita")

cap_z <- "Multilevel models of pro-immigrant attitudes. Random intercepts at the country level. Higher values of the DV indicate more pro-immigrant attitudes."



#Out to Latex
texreg(list(m1z, m2z, m3z, m4z),
       custom.coef.map = coef_labs_z,
       stars = c(0.1, 0.05, 0.01),
       caption = cap_z,
       file = "tables/z_score_models.tex",
       label = "tab_z_models")


```


# Models with Additional Country-Level Vars 

```{r estimates-with-controls}

# Interactions with Log GDP/Capita
gdp_cosmo <- lmer(imm_index ~ close_world * refugee_change +  ln_gdppc *  close_world + edu_cat + trust_generalized + 
             income_log + urban + unemp + ln_age +  female + unemp +  (1| iso3c), data = wvs_latam) 

gdp_nat <- lmer(imm_index ~ natpride * refugee_change + ln_gdppc * natpride + edu_cat + trust_generalized + income_log +
                    urban + ln_age + female + unemp + age + female + (1| iso3c), data = wvs_latam)


# Interactions with Unemployment 
unemp_cosmo <- lmer(imm_index ~ close_world * refugee_change + unemploytotal *  close_world + edu_cat + trust_generalized + 
             income_log + urban + unemp + ln_age +  female + unemp +  (1| iso3c), data = wvs_latam) 

unemp_nat <- lmer(imm_index ~ natpride * refugee_change + unemploytotal * natpride + edu_cat + trust_generalized + income_log +
                    urban + ln_age + female + unemp + age + female + (1| iso3c), data = wvs_latam)


# Interactions with Polyarchy 
polyarchy_cosmo <- lmer(imm_index ~ close_world * refugee_change +  v2x_polyarchy *  close_world + edu_cat + trust_generalized + 
             income_log + urban + unemp + ln_age +  female + unemp +  (1| iso3c), data = wvs_latam) 

polyarchy_nat <- lmer(imm_index ~ natpride * refugee_change + v2x_polyarchy * natpride + edu_cat + trust_generalized + income_log +
                    urban + ln_age + female + unemp + age + female + (1| iso3c), data = wvs_latam)


coef_labs_alt <- list(
                  "close_world" = "Cosmopolitanism", 
                   "close_world:refugee_change" = "a $\\times$ Refugee Change",
                   "close_world:ln_gdppc" = "b $\\times$ Log GDP/Capita",
                   "close_world:unemploytotal" = "c $\\times$ Unemployment",
                   "close_world:v2x_polyarchy" = "g $\\times$ Polyarchy",
                  "natpride" = "National Pride",
                  "natpride:refugee_change" = "d $\\times$ Refugee Change",
                  "natpride:ln_gdppc" = "e $\\times$ Log GDP/Capita",
                  "natpride:unemploytotal" = "f $\\times$ Unemployment",
                  "natpride:v2x_polyarchy" = "h $\\times$ Polyarchy",
                  "trust_generalized" = "Interpersonal Trust",
                  "edu_catSecondary" = "Secondary Edu", 
                  "edu_catPost-Secondary" = "Post-Secondary Edu", 
                  "ln_age" = "Age (logged)",
                  "urban1" = "Urban",
                  "unemp1" = "Unemployed",
                  "female1" = "Female",
                  "income_log" = "Individual Income", 
                  "refugee_change" = "Refugee Change", 
                   "ln_gdppc" = "Log GDP/Capita",
                  "unemploytotal" = "Unemployment",
                  "v2x_polyarchy" = "Polyarchy")


texreg(list(gdp_cosmo, unemp_cosmo, polyarchy_cosmo),  
       stars = c(0.1, 0.05, 0.01), 
       custom.coef.map = coef_labs_alt,
       file = "tables/si_control_cosmo.tex",
       label = "tab:si_control_mods")


texreg(list(gdp_nat, unemp_nat, polyarchy_nat),  
       stars = c(0.1, 0.05, 0.01), 
       custom.coef.map = coef_labs_alt,
       file = "tables/si_control_nat.tex",
       label = "tab:si_control_mods")



```


# Models by Urban/Rural

```{r urban-rural-mods}

# Make urban a factor for plotting
wvs_latam <- wvs_latam %>% 
  mutate(urban_fac = case_match(urban, 
                               "1" ~ "Urban", 
                               "0" ~ "Rural"))

# Conditional Effect of Cosmopolitanism 
m2_ur <- lmer(imm_index ~ close_world * refugee_change * urban_fac + edu_cat + trust_generalized + 
             income_log + unemp + ln_age +  female + unemp + ln_gdppc + (1| iso3c), data = wvs_latam) 

# plot
ur1 <- plot_slopes(m2_ur, variables = "close_world", condition = c("refugee_change", "urban_fac")) +
  facet_grid(~urban_fac) +
  geom_hline(yintercept = 0, linetype = "dashed", color = "red") +
  labs(x = "% Change in Refugee Population", y = "Cosmopolitanism") +
  theme_pubr(base_size = 13) +
  theme(legend.position = "none",  strip.text = element_text(face = "bold"))


# Conditional Effect of National Pride
m4_ur <- lmer(imm_index ~ natpride * refugee_change * urban_fac +  edu_cat + trust_generalized + income_log
            + ln_age + female + unemp + ln_gdppc + age + female + (1| iso3c), data = wvs_latam)


ur2 <- plot_slopes(m4_ur, variables = "natpride", condition = c("refugee_change", "urban_fac")) +
  facet_grid(~urban_fac) +
  geom_hline(yintercept = 0, linetype = "dashed", color = "red") +
  labs(x = "% Change in Refugee Population", y = "Nationalism") +
  theme_pubr(base_size = 13) +
  theme(legend.position = "none", strip.text = element_text(face = "bold"))

# combine plots
ggarrange(ur1 + rremove("x.title"), ur2, 
          ncol = 1, align = "v") |>  
  annotate_figure(top = text_grob("Effect on Pro-Immigration Attitudes",
                                  face = "bold", size = 14))

# save
ggsave("figures/ur_plot.pdf", height = 6, width = 8)


```


# Models by News Consumption

```{r news-models}

# Cosmo model
news_cosmo <- lmer(imm_index ~ close_world * refugee_change * news_consumption + edu_cat + trust_generalized + 
                      income_log + unemp + ln_age +  female + unemp + ln_gdppc + (1| iso3c), data = wvs_latam) 

# Summary stats in sample
summary(news_cosmo@frame$news_consumption)

# Get slopes
news_slopes_c <-  avg_slopes(news_cosmo, 
                             variables = "close_world", 
                             by = c("refugee_change", "news_consumption")) 

# Label for plotting
news_slopes_c <-  news_slopes_c |> 
  mutate(
    news = case_match(news_consumption,
                      0    ~ "Min. News",
                      2.5  ~ "Median News", 
                      4    ~ "Max. News"),
    news = factor(news, levels = c("Min. News", "Median News", "Max. News"))
  ) 

# Plotting
plot_c <- ggplot(na.omit(news_slopes_c), aes(x = refugee_change, y = estimate, fill = news)) +
  geom_line() +
  geom_ribbon(aes(ymin = conf.low, ymax = conf.high), alpha = .2) +
  facet_grid(~news) + 
  geom_hline(yintercept = 0, linetype = "dashed", color = "red") +
  labs(y = "Effect of Cosmopolitanism") +
  theme_pubclean() +
  theme(legend.position = "none") 

# Nationalism model
news_nat <- lmer(imm_index ~ natpride * refugee_change * news_consumption + edu + trust_generalized + 
                    income_log + unemp + ln_age +  female + unemp + ln_gdppc + (1| iso3c), data = wvs_latam)

# Get slopes 
summary(news_nat@frame$news_consumption)

# Get slopes
news_slopes_n <-  avg_slopes(news_nat, 
                             variables = "natpride", 
                             by = c("refugee_change", "news_consumption")) 

news_slopes_n <-  news_slopes_n |> 
  
  mutate(
    news = case_match(news_consumption,
                      0.000    ~ "Min. News",
                      2.375 ~ "Median News", 
                      4.000    ~ "Max. News"),
    news = factor(news, levels = c("Min. News", "Median News", "Max. News"))
  ) 
  

# Plotting
plot_n <-  ggplot(na.omit(news_slopes_n), aes(x = refugee_change, y = estimate, fill = news)) +
  geom_line() +
  geom_ribbon(aes(ymin = conf.low, ymax = conf.high), alpha = .2) +
  facet_grid(~news) + 
  geom_hline(yintercept = 0, linetype = "dashed", color = "red") +
  labs(y = "Effect of Nationalism", x = "% Change in Refugee Population") +
  theme_pubclean() +
  theme(legend.position = "none") 


ggarrange(plot_c + rremove("x.title"), 
          plot_n + rremove("x.title"), 
          ncol = 1, align = "v") |>  
  annotate_figure(top = text_grob("Effect on Pro-Immigration Attitudes",
                                  face = "bold", size = 14))


ggsave("figures/cond_news_plot.pdf", height = 6, width = 8)



```


# Models with Percent Change in Migration  

```{r estimate-mig-mods}

# Main Effect of Cosmo
m1_mig <- lmer(imm_index ~ close_world + edu_cat + trust_generalized + income_log + urban +
             ln_age + female + unemp + mig_change + ln_gdppc +  (1| iso3c), data = wvs_latam) 


# Conditional Effect of Cosmopolitanism 
m2_mig <- lmer(imm_index ~ close_world * mig_change + edu_cat + trust_generalized + 
             income_log + urban + unemp + ln_age +  female + unemp + ln_gdppc +  (1| iso3c), data = wvs_latam) 


# Main Effect of National Pride
m3_mig <- lmer(imm_index ~ natpride + edu_cat + trust_generalized + income_log + ln_age + female + 
             unemp + urban + mig_change + ln_gdppc + (1| iso3c), data = wvs_latam) 


# Conditional Effect of National Pride
m4_mig <- lmer(imm_index ~ natpride * mig_change + edu_cat + trust_generalized + income_log + urban +
             ln_age + female + unemp + ln_gdppc + age + female + (1| iso3c), data = wvs_latam)

# Print the results
coef_labs_mig <- list(
                  "close_world" = "Cosmopolitanism", 
                   "close_world:mig_change" = "Cosmopolitanism $\\times$ Migrant Change",
                  "natpride" = "National Pride",
                  "natpride:mig_change" = "National Pride $\\times$ Migrant Change",
                  "trust_generalized" = "Interpersonal Trust",
                  "edu_catSecondary" = "Secondary Edu", 
                  "edu_catPost-Secondary" = "Post-Secondary Edu", 
                  "ln_age" = "Age (logged)",
                  "urban1" = "Urban",
                  "unemp1" = "Unemployed",
                  "female1" = "Female",
                  "income_log" = "Individual Income", 
                  "mig_change" = "Migrant Change", 
                   "ln_gdppc" = "Log GDP/Capita")


texreg(list(m1_mig, m2_mig, m3_mig, m4_mig),
       custom.coef.map = coef_labs_mig,
       stars = c(0.1, 0.05, 0.01),
       file = "tables/migrant_change_models.tex")



```


# Predict Immig with Close Measures


```{r construct}

# Select immigration index and the "close to" variables 
wvs_close <- wvs_latam %>% 
  select(imm_index, close_town, close_region, close_country, close_continent, close_world)

# Regress imm Index on all Variables
a1 <- lm(imm_index ~ ., data = wvs_close)

# Create a dot and whisker plot
dw_plot_a1 <- dwplot(a1)

# Make some labels 
close_labs <- c("Close to World", "Close to Continent",
                "Close to Country", "Close to Region", "Close to Town") 

# Plot
dw_plot_a1 +
  scale_y_discrete(labels = close_labs) +
  labs(title = "DV: Pro-Immigration Index", x = "Coefficient Estimate") +
  geom_vline(xintercept = 0, linetype = "dashed") +
  theme_pubr(base_size = 14) +
  theme(axis.title = element_text(face = "bold"), 
        legend.position = "none")

ggsave(file = here("figures/cosmo_dot_plot.pdf"),
       height = 6,
       width = 8)

```


# Analysis with Full WVS Sample 

```{r load-data}

# Read in the WVS data for the full sample
wvs_full <- read_csv("wvs_all.csv")


```


### Recoding/Merging for Models 

```{r shaping-full}
# Make Education a Factor
wvs_full <- wvs_full %>% 
  mutate(edu_cat = case_match(edu,
                              1:2 ~ "Primary", 
                              1:4 ~ "Secondary",
                              5:8 ~ "Post-Secondary"))

wvs_full$edu_cat <- factor(wvs_full$edu_cat, levels = c("Primary", 
                                              "Secondary", 
                                              "Post-Secondary"))


# Log the Income Level Variable 
wvs_full$income_log <- log(wvs_full$income_lvl)

# Log Age 
wvs_full$ln_age <- log(wvs_full$age)

# Make a couple of variables factors
wvs_full$unemp <- factor(wvs_full$unemp) #unemployment
wvs_full$female <- factor(wvs_full$female) #female
wvs_full$urban <-  factor(wvs_full$urban) # urban rural


```

## Models Full Sample 

- Multilevel models with country-level random effects. 
- First estimate main effect of cosmopolitanism and national pride, then estimate interaction with change in number of refugees. 


```{r estimate-full-sample}

# Main Effect of Cosmo
m1_full <- lmer(imm_index ~ close_world + edu_cat + trust_generalized + income_log + urban +
             ln_age + female + unemp + refugee_change + ln_gdppc +  (1| iso3c), data = wvs_full) 


# Conditional Effect of Cosmopolitanism 
m2_full <- lmer(imm_index ~ close_world * refugee_change + edu_cat + trust_generalized + 
             income_log + urban + unemp + ln_age +  female + unemp + ln_gdppc +  (1| iso3c), data = wvs_full) 


# Main Effect of National Pride
m3_full <- lmer(imm_index ~ natpride + edu_cat + trust_generalized + income_log + ln_age + female + 
             unemp + urban + refugee_change + ln_gdppc + (1| iso3c), data = wvs_full) 



# Conditional Effect of National Pride
m4_full <- lmer(imm_index ~ natpride * refugee_change + edu_cat + trust_generalized + income_log + urban +
             ln_age + female + unemp + ln_gdppc + age + female + (1| iso3c), data = wvs_full)


# List of our models 
full_mods <- list(m1_full, m2_full, m3_full, m4_full)

# Coefficient Labels 
coef_labs <- list(
                  "close_world" = "Cosmopolitanism", 
                   "close_world:refugee_change" = "Cosmopolitanism $\\times$ Refugee Change",
                  "natpride" = "National Pride",
                  "natpride:refugee_change" = "National Pride $\\times$ Refugee Change",
                  "trust_generalized" = "Interpersonal Trust",
                  "edu_catSecondary" = "Secondary Edu", 
                  "edu_catPost-Secondary" = "Post-Secondary Edu", 
                  "ln_age" = "Age (logged)",
                  "urban1" = "Urban",
                  "unemp1" = "Unemployed",
                  "female1" = "Female",
                  "income_log" = "Individual Income", 
                  "refugee_change" = "Refugee Change", 
                   "ln_gdppc" = "Log GDP/Capita")

# Caption 
cap <- "Multilevel models of pro-immigrant attitudes. Random intercepts at the country level. Higher values of the DV indicate more pro-immigrant attitudes."

#Out to Latex
texreg(full_mods,
       custom.coef.map = coef_labs,
       stars = c(0.1, 0.05, 0.01),
       caption = cap,
       file = "tables/world_sample_mods.tex",
       label = "tab:world_sample")



```


# Model Each Item in Immigration Index Seperately

```{r break-up-index}

# Grab the DVS
dvs <- wvs_latam |> 
  select(starts_with("imm")) |> 
  colnames()

# Six items in the index
dvs <- dvs[1:6]

# Labels for each
labels <- c(
  "Development",
  "Diversity",
  "Crime Rate",
  "Risk of Terrorism",
  "Unemployment",
  "Social Conflict"
)

# save results to lists
cosmo_results <- list()
nat_results <- list()


# Iterate over dependent variables
for (i in 1:length(dvs)) { # begin loop
  
  # Define formulas
  cosmo_formula <- as.formula(paste(dvs[i], "~ close_world * refugee_change + edu_cat + trust_generalized + income_log + urban + unemp + ln_age + female + unemp + ln_gdppc + (1|iso3c)"))
  nat_formula <- as.formula(paste(dvs[i], "~ natpride * refugee_change + edu_cat + trust_generalized + income_log + urban + unemp + ln_age + female + unemp + ln_gdppc + (1|iso3c)"))
  
  # Fit models
  cosmo_fit <- lmer(cosmo_formula, data = wvs_latam) # Cosmopolitanism model
  nat_fit <- lmer(nat_formula, data = wvs_latam) # Nationalism model
  
  # Plot
  cosmo_plot <- plot_slopes(cosmo_fit, variables = "close_world", condition = "refugee_change") +
    labs(title = labels[i]) + 
    geom_hline(yintercept = 0, linetype = "dashed", color = "red") + 
    scale_y_continuous(limits = c(-.07, 0.08)) +
    theme_pubr(base_size = 11) +
    theme(axis.title = element_blank())
  
  cosmo_plot <- plot(cosmo_plot)
  
  nat_plot <- plot_slopes(nat_fit, variables = "natpride", condition = "refugee_change") +
    labs(title = labels[i]) +
    geom_hline(yintercept = 0, linetype = "dashed", color = "red") + 
    scale_y_continuous(limits = c(-.13, 0.08)) + 
    theme_pubr(base_size = 11) +
    theme(axis.title = element_blank())
  
  nat_plot <- plot(nat_plot)
  
  # Save results
  cosmo_results[[i]] <- cosmo_plot 
  nat_results[[i]] <- nat_plot 
  
} # end loop

# combine cosmo plots
ggarrange(cosmo_results[[1]],
          cosmo_results[[2]] + rremove("y.text"), 
          cosmo_results[[3]] + rremove("y.text"),
          cosmo_results[[4]],
          cosmo_results[[5]] + rremove("y.text"), 
          cosmo_results[[6]]+ rremove("y.text")) |>  
  # annotate
  annotate_figure(
    left = text_grob("Effect of Cosmopolitanism", rot = 90, size = 14, face = "bold"),
    bottom = text_grob("% Change in Refugee Population", size = 14, face = "bold")
    )

# save
ggsave("figures/cosmo_disagg.pdf", 
       height = 6,
       width = 8)

# combine nat plots
ggarrange(
  nat_results[[1]],
  nat_results[[2]] + rremove("y.text"), 
  nat_results[[3]] + rremove("y.text"),
  nat_results[[4]], 
  nat_results[[5]] + rremove("y.text"), 
  nat_results[[6]] + rremove("y.text")) |>  
  # annotate 
    annotate_figure(
      left = text_grob("Effect of Nationalism", rot = 90, size = 14, face = "bold"),
      bottom = text_grob("% Change in Refugee Population", size = 14, face = "bold")
      )

# save
ggsave("figures/nat_disagg.pdf", 
       height = 6,
       width = 8)


```


# Control for and Prediction Nationalism/Cosmopolitanism

```{r mediate}

# predict cosmopolitanism
imm_cosmo <- lmer(close_world ~  + edu_cat + trust_generalized + income_log + urban +
             ln_age + female + unemp + refugee_change + ln_gdppc +  (1| iso3c), data = wvs_latam)

# predict nationalism 
imm_nat <- lmer(natpride ~ edu_cat + trust_generalized + income_log + ln_age + female
                 + unemp + urban + refugee_change + ln_gdppc + (1| iso3c), data = wvs_latam)


# control for mean cosmo
control_cosmo <- lmer(imm_index ~ close_world * refugee_change + edu_cat + trust_generalized + 
             income_log + urban + unemp + ln_age +  female + unemp + 
               ln_gdppc + mean_cosmo + (1| iso3c), data = wvs_latam) 

# control for mean nationalism
control_nat <- lmer(imm_index ~ natpride * refugee_change + edu_cat + trust_generalized + income_log + ln_age + female + 
             unemp + urban + ln_gdppc + mean_pride + (1| iso3c), data = wvs_latam) 

these_mods <- list(imm_cosmo, imm_nat, control_cosmo, control_nat)

these_coef_labs <- list(
                  "close_world" = "Cosmopolitanism", 
                   "close_world:refugee_change" = "Cosmopolitanism $\\times$ Refugee Change",
                  "natpride" = "National Pride",
                  "natpride:refugee_change" = "National Pride $\\times$ Refugee Change",
                  "trust_generalized" = "Interpersonal Trust",
                  "edu_catSecondary" = "Secondary Edu", 
                  "edu_catPost-Secondary" = "Post-Secondary Edu", 
                  "ln_age" = "Age (logged)",
                  "urban1" = "Urban",
                  "unemp1" = "Unemployed",
                  "female1" = "Female",
                  "income_log" = "Individual Income", 
                  "refugee_change" = "Refugee Change", 
                  "ln_gdppc" = "Log GDP/Capita", 
                  "mean_pride" = "National Pride (country mean)",
                  "mean_cosmo" = "Cosmopolitanism (country mean)")



texreg(these_mods, 
       custom.coef.map = these_coef_labs, 
       stars = c(0.1, 0.05, 0.01),
       file = "tables/country_control_mods.tex")


```




















